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HIGHLIGHTS 


•  The  experiment  and  simulation  of  LIB  nail  penetration  test  are  performed. 

•  We  model  the  behavior  of  short-circuit  LIBs  using  porous  electrode  theory. 

•  Accurately  describing  the  short-circuit  voltage  is  essential. 

•  Nail  penetration  induced  thermal  runaway  is  described  by  abuse  equations. 

•  The  short-circuit  heating  determines  the  moment  of  thermal  runaway. 
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Nail  penetration  into  a  battery  pack,  resulting  in  a  state  of  short-circuit  and  thus  burning,  is  likely  to 
occur  in  electric  car  collisions.  To  demonstrate  the  behavior  of  a  specific  battery  when  subject  to  such 
incidents,  a  standard  nail  penetration  test  is  usually  performed;  however,  conducting  such  an  experi¬ 
ment  is  money  consuming.  The  purpose  of  this  study  is  to  propose  a  numerical  electrochemical  model 
that  can  simulate  the  test  accurately.  This  simulation  makes  two  accurate  predictions.  First,  we  are  able  to 
model  short-circuited  lithium-ion  batteries  (LIBs)  via  electrochemical  governing  equations  so  that  the 
mass  and  charge  transfer  effect  could  be  considered.  Second,  the  temperature  variation  of  the  cell  during 
and  after  nail  penetration  is  accurately  predicted  with  the  help  of  simulating  the  temperature  distri¬ 
bution  of  thermal  runaway  cells  by  thermal  abuse  equations.  According  to  this  nail  penetration  model, 
both  the  onset  of  battery  thermal  runaway  and  the  cell  temperature  profile  of  the  test  are  obtained,  both 
of  which  are  well  fitted  with  our  experimental  results. 

©  2013  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

Lithium-ion  batteries  (LIBs)  are  widely  used  in  portable  elec¬ 
tronic  devices.  With  characteristics  of  high  power  density  and  little 
memorizing  effect,  LIBs  are  now  also  being  utilized  as  the  power 
source  for  electric-powered  transportation  tools,  such  as  EVs 
(electric  vehicles)  and  HEVs  (hybrid  electric  vehicles)  [1,2].  Since 
passenger  safety  is  closely  related  to  the  vehicle’s  design,  a  thor¬ 
ough  understanding  of  the  EV’s  and  HEV’s  LIB  behavior  under 
extreme  environment  is  essential,  especially  in  abuse  conditions 
which  can  potentially  cause  inflammation  or  explosion  of  the  bat¬ 
tery  pack.  As  a  result,  various  tests  can  be  conducted  to  give  an 
insight  of  the  battery’s  property  under  different  abuse  conditions. 
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In  this  paper,  we  focus  on  one  test— the  nail  penetration  test 
[3-6  —for  battery  safety.  By  driving  a  steel  nail  directly  into  the 
LIB,  the  negative  and  positive  electrodes  are  partially  connected 
inside  the  jelly-roll  by  the  electric  conducting  nail.  At  the  point  of 
penetration,  electrons  can  flow  freely  between  the  two  elec¬ 
trodes  through  the  nail,  resulting  in  a  local  internal  short-circuit 
that  occurs  at  the  nail  and  jelly-roll  interface.  Joule  and  entropic 
heat  are  thus  generated  due  to  the  large  electric  current  density 
and  the  redox  reaction  of  lithium,  respectively  [7].  This  excessive 
heat  will  propagate  throughout  the  cell  from  the  point  of  pene¬ 
tration,  raising  the  temperature  of  the  jelly-roll  within.  If  no 
current  shut-down  or  temperature  suppressing  device  is  effec¬ 
tive,  the  cell  temperature  could  upsurge  high  enough  to  trigger 
the  thermal  decomposing  reactions  of  the  jelly-roll  materials. 
Once  the  reactions  are  activated,  it  is  generally  inescapable  for 
battery  thermal  runaway  to  occur  8]. 

By  conducting  the  battery  nail  penetration  test,  valuable  infor¬ 
mation  about  the  behavior  of  batteries  under  internal  short-circuit 
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could  be  provided.  However,  such  experiments  may  be  money 
consuming,  and  with  the  potential  of  battery  explosion,  nail 
penetration  must  take  place  inside  a  fire-proof  infrastructure, 
causing  it  a  dangerous  test  to  conduct. 

Besides  experiments,  numerical  methods  provide  a  sensible 
alternative  means  to  save  money.  Modeling  the  nail  penetration  in¬ 
cludes  the  simulation  of  (i)  LIB  short-circuit  and  (ii)  material  thermal 
decomposition.  When  encountering  LIB  short-circuit,  the  heat  gen¬ 
eration  rate  can  be  simulated  by  considering  only  the  ohmic  heating 
of  the  cell  due  to  a  large  current  density  discharge  9,10].  However,  a 
more  sophisticated  method  is  to  consider  the  entropic  heating  as 
well,  and  the  temperature  of  the  battery  can  be  simulated  if  the 
thermal  properties  of  the  cell  are  known  11  -14].  On  the  other  hand,  a 
crush  experiment  could  be  performed  to  understand  the  effect  of 
internal  short-circuit  due  to  physical  abuse  15].  The  experiment  re¬ 
veals  the  situation  when  LIBs  are  dropped  or  bent  by  an  outer  force. 
Another  condition  of  LIB  short-circuit  is  the  contamination  of  metal 
particle  during  manufacturing,  and  studies  have  shown  that  there  is  a 
significant  difference  between  the  cell  behaviors  of  nail  penetration 
induced  and  contamination  induced  short-circuit  16]. 

In  this  work,  Section  2  presents  the  setup  and  the  experimental 
results  of  the  nail  penetration  test.  Section  3  describes  the 
improvement  of  battery  short-circuit  simulation  by  simulating  the 
battery  short-circuit  process  via  electrochemical  governing  equa¬ 
tions  and  by  considering  both  the  Joule  and  entropic  heat  generation 
rates.  In  so  doing,  the  initial  heating  of  the  battery  after  penetration, 
which  accounts  for  the  triggering  of  thermal  runaway,  can  be 
modeled  by  a  more  rigorous  treatment.  Moreover,  the  influence  of 
cell  parameters  on  battery  heating,  such  as  the  state  of  charge  (SOC), 
could  also  be  qualitatively  and  quantitatively  analyzed.  In  addition 
to  battery  short-circuit,  we  adopted  thermal  abuse  equations  [  17,18  ], 
described  in  Section  4,  to  accurately  model  the  exothermic  material 
decomposition  of  LIB  materials  when  exposed  to  high  temperature. 
We  believe  an  exact  characterization  of  the  reaction’s  heat  genera¬ 
tion  rate  is  essential  for  evaluating  the  effectiveness  of  a  specific 
battery  pack  nail  penetration  safety  design.  In  Section  5,  we  combine 
the  LIB  short-circuit  simulation  approach  of  Section  3  and  the 
thermal  abuse  equations  of  Section  4  to  construct  the  LIB  nail 
penetration  model.  Section  6  is  a  conclusion. 

2.  Nail  penetration  test  of  LIBs 

The  experiment  setup  of  the  nail  penetration  test  is  shown  in 
Fig.  1.  The  type  of  LIB  we  used  in  the  test  is  an  ITRI  (Industrial 
Technology  Research  Institute  of  Taiwan)  5.25  Ah  lithium  manga¬ 
nese  oxide  (LMO)  soft  package  prismatic  cell,  initially  charged  to 


100%  SOC  before  penetration.  The  cell  was  fresh,  and  no  additional 
safety  designs  to  prevent  the  battery  from  short-circuit  or  thermal 
runaway  were  embedded  within.  A  total  of  five  thermal  couples, 
mounted  at  locations  as  indicated  in  Fig.  1(a),  were  used  to  record 
the  battery  surface  temperature.  The  picture  of  the  actual  battery 
and  thermal  couple  layout  is  shown  in  Fig.  1(b).  A  3  mm  diameter 
stainless  steel  nail  was  used  as  the  penetrating  nail.  The  point  of 
penetration  was  located  at  the  center  of  the  battery,  with  the 
penetrating  direction  being  perpendicular  to  the  battery  surface. 
The  nail  was  programmed  to  penetrate  through  the  cell  with  the 
speed  of  10  mm  s-1,  and  remains  in  the  cell  throughout  the 
experiment.  The  whole  nail  penetration  test  was  conducted  inside 
an  explosion-proof  chamber,  with  ventilation  turned  on  after 
penetration  occurred.  Initial  temperature  of  the  battery  and  the 
ambient  was  kept  at  25.7  °C.  Due  to  a  possibility  of  battery  dilation, 
two  fixed  Bakelites  were  in  direct  contact  with  the  top  and  bottom 
surfaces  of  the  cell  to  prevent  the  thermal  couples  from  detaching 
during  the  test. 

The  temperature  evolution  recorded  by  each  thermal  couple  is 
shown  in  Fig.  2,  where  there  is  a  sudden  jump  of  cell  temperature  to 
over  300  °C  at  approximately  50  s  after  nail  penetration.  This 
abrupt  rise  of  temperature  indicates  the  occurrence  of  battery 
thermal  runaway.  Furthermore,  thermal  couple  E  captured  the 
greatest  temperature  peak  value  of  the  five,  while  the  two  smallest 
peaks  were  recorded  by  thermal  couples  C  and  D.  Near  the  center  of 
the  cell,  where  short-circuit  took  place,  thermal  couple  E  received 
the  maximum  heat  from  short-circuit  and  minimum  cooling  from 
ventilation.  In  contrast,  the  peaks  were  smallest  in  thermal  couples 
C  and  D  as  they  were  located  far  away  from  the  point  of  penetration 
and  near  the  edge  of  the  battery,  where  convection  had  the 
strongest  effect. 

3.  LIB  short-circuit  simulation 

For  simulating  LIBs  operating  under  normal  conditions,  the 
porous  electrode  theory  19-24],  presented  in  the  following 
paragraph  with  a  concise  description,  provides  a  background  of  the 
battery’s  electric  and  mass  transfer  process.  The  theory  usually 
treats  the  temperature  as  a  constant  parameter  and  simplifies  the 
jelly-roll  of  the  cell  to  a  one-dimensional  geometry  that  is 
composed  of  five  segments,  which  correspond  to  the  thickness  of 
the  two  current  collectors  and  electrodes  and  the  separator.  Fig.  3  is 
an  example  of  the  structure,  where  L  is  the  cell  thickness  and  the 
leftmost  boundary  of  the  negative  electrode  current  collector  is  set 
as  the  origin  of  the  x*-axis,  which  is  the  theory’s  one-dimensional 
coordinate. 


99  mm 


Fig.  1.  The  layout  of  the  ITRI  5.25  Ah  soft  package  prismatic  cell  nail  penetration  test  with  the  upper  Bakelite  removed  for  a  clearer  observation,  (a)  The  cell  has  a  99  x  130  mm2 
cross-sectional  area  and  5  mm  thick.  The  mounting  coordinates,  with  the  origin  at  the  center  of  the  cell,  for  thermal  couples  A,  B,  C,  D,  and  E  are  (-2,  2),  (-2,  -2),  (3.5,  5),  (3.5,  -5), 
and  (1,  0),  respectively,  (b)  A  photo  image  of  the  actual  test  setup  with  the  five  mounted  thermal  couples. 
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Fig.  2.  Experimental  result  of  the  ITRI  5.25  Ah  soft  package  prismatic  cell  nail  pene¬ 
tration  test.  The  five  temperature  curves  are  obtained  by  the  measurement  of  the 
corresponding  thermal  couples  shown  in  Fig.  1(a).  The  inset  is  a  photograph  of  the  test 
setup  after  the  penetration.  Although  the  cell  generated  large  internal  pressure  due  to 
thermal  runaway  and  broke  the  upper  Bakelite,  the  five  thermal  couples  remained 
attached  to  the  cell  properly. 


the  Butler-Volmer  equation,  which  describes  the  difference  be¬ 
tween  reduction  and  oxidation  reaction  rates  of  lithium  at  the 
electrode-electrolyte  interface: 


J\A,i 


Cf)°'5(Csf)a5 


x 


FVs  ,i 

2  RT 


-  exp 


-FVs  ,i 

2  RT 


i  =  n,p, 


(4) 


where  r\s  i  =  cf>s  i  -  cf>{  i  -  Uf  -  FjUiRfi  is  the  surface  over-potential. 
The  mass  transfer  process  inside  electrolyte  phase  in  the  sepa¬ 
rator  region,  negative  electrode,  and  positive  electrode  includes  the 
diffusion  of  lithium-salt,  cations,  and  anions: 


dCi 

^at 


Dffv*q 


/MV*t+ 


fli/LU 


1  -t 


i  =  e,n,p, 


(5) 


where  lithium-ion  diffusivity  of  the  electrolyte  phase  is  corrected 
by  the  Bruggeman  exponential  y  as  Dff  =  DyeJ^  Potential  distri¬ 
bution  across  the  electrolyte  phase  with  respect  to  its  current 
density  is  evaluated  by  the  modified  Ohm’s  law: 


From  the  porous  electrode  theory,  Fields  second  law  of  diffusion 
is  applied  to  evaluate  the  change  of  lithium  concentration  due  to 
solid  phase  diffusion  inside  the  negative  and  positive  electrodes  of 
the  battery  [25-27  : 


at 


^  csj  2  0cs  j 
dr 2  +  r  6r 


i  =  n,p, 
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2Kf?RT  f  x 

lu  =  -*uV  *.*+ — p — (l  -t+)v  Inc,  i  =  e,  n,p.  (6) 

An  alternative  form  of  Faraday’s  law,  similar  to  Equation  (3),  is 
used  to  link  the  electrolyte  phase  current  density  gradient  and  the 
pore  wall  lithium-ion  flux: 

V*/,j  =  FcyLU,  i  =  n,p.  (7) 


where  r  is  the  radial  coordinate  of  the  electrode’s  particle.  Ohm’s 
law  describes  the  connection  between  solid  phase  current  density 
and  the  solid  phase  electric  potential  gradient: 

4,i  =  V*4,,  i  =  n,p,  (2) 

where  the  effective  solid  phase  electric  conductivity  is  calculated 
by  afj  =  <7S  j(l  -  eSji  -  £f  j),  and  V*  =  d/dx*  is  the  gradient  operator 
with  respect  to  the  model’s  one-dimensional  coordinate.  An  alter¬ 
native  form  of  Faraday’s  law  relates  the  solid  phase  current  density 
gradient  and  the  pore  wall  lithium-ion  flux: 


V*/s i  =  -Fajn  i  =  n,p,  (3) 

where  F  is  the  Fdraday  constant  and  a*  is  the  specific  interfacial  area, 

determined  by  az-  =  3(1  -  eSti  -  ei ,i)lRs,i  with  the  radius  of  electrode 

particle  Rs.  The  magnitude  and  direction  of  the  flux  is  calculated  by 


3.1.  Short-circuit  simulation  setup 

Short-circuit  of  battery  indicates  that  the  positive  electrode’s 
electric  potential  value  (0S,P)  is  very  close  to  the  negative  elec¬ 
trode’s  (0 s,n).  To  simulate  short-circuit,  we  modified  the  porous 
electrode  theory  by  artificially  giving  a  potential  value  at  the 
positive  end  of  the  battery  (x*  =  L),  which  serves  as  a  boundary 
condition  in  the  numerical  model,  i.e.,  </>s,p(x*5  01  -x*=l  =  0s,p(OO 
with  0S)P (L,t)  being  given.  Initially,  for  the  simplest  case  possible, 
we  assumed  that  </>SiP(L,t)  =  0.5  V  (relative  to  ground)  and  it 
would  be  kept  constant  throughout  the  whole  short-circuit 
process.  However,  previous  experiment  results  showed  that  0S,P 
decayed  with  time  soon  after  short-circuit  occurred  [28,29]. 
Therefore,  we  reconsidered  0sp(L,t)  as  0s>p(L,t)  =  0.5exp 
(-0.07 1)  V,  which  is  a  close  fit  of  electric  potential  to  the 
experimentally  measured  time  decaying  profile. 


Negative  Porous  Negative 
Current  Collector  Electrode 


Separator 


Porous  Positive  Positive 
Electrode  Current  Collector 


Fig.  3.  The  one-dimensional  geometry  in  the  porous  electrode  theory. 
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Battery  short-circuit  causes  a  severe  rise  of  cell  temperature.  To 
simulate  the  corresponding  heating  of  the  battery,  we  utilized  the 
simplified  form  of  LIB  heat  generation  rate  provided  by  Bernardi 
et  al.  [30]: 


S 


sc 


l 


^  aih,i(Uf  0s, i) 


^ih,i 


i  =  n,p,  (8) 


which  was  frequently  adopted  in  LIB  thermal  analysis  [8,31-37]. 
There  are  two  major  heat  sources  produced  by  the  electrochemical 
process  of  LIBs:  the  irreversible  Joule  heating  and  the  reversible 
entropic  heating.  They  are  represented  by  the  first  and  second  term 
of  Equation  (8),  respectively.  The  entropy  of  the  carbon  negative 
electrode  Li*C6  is  provided  by  Srinivasan  et  al.  38]  in  the  form  of: 


0s, n*  0s, p>  01, n>  01, p>  01, s>  Q,n>  Q,p>  On  Cp,  cs,  jm  Jp*  f^n*  ^d  Up)  in  the 
porous  electrode  theory  with  the  help  of  numerical  iteration.  The 
cell  capacity  is  obtained  from  the  product  of  the  applied  discharge 
current  density  Is  (A  nrT2),  total  effective  interfacial  area  of  cell  Aceii 
(m2),  and  discharge  time  (in  hours),  where  the  1C  discharge  current 
density,  9.77  Am-2,  is  the  value  of  Is  when  the  vertical  drop  of  the 
simulated  voltage  occurs  at  precisely  3600  s.  The  value  of  Aceii  could 
then  be  determined  from  dividing  the  cell’s  nominal  capacity 
(5.25  Ah)  by  the  1C  discharge  current  density  and  1  h,  gaining  the 
value  of  Aceii  =  0.537  m2.  Our  simulation  result  is  shown  Fig.  4,  and 
the  corresponding  parameters  and  electrodes’  open-circuit  poten¬ 
tials  (OCPs)  are  respectively  presented  in  fable  1  and  in  Fig.  5. 


3.3.  Simulation  result 


d [/<-  _  344.1347148  exp(-32.9633287x  + 8.316711484) 

dT^  ~  1  +  749.0756003  exp(-34.79099646x  +  8.887143624) 

-  0.8520278805*  +  0.362299229*2  +  0.2698001697, 

(9) 

and  the  entropy  of  the  LiyMn204  positive  electrode  is  given  as: 

It  should  be  noted  that  the  unit  of  the  two  electrodes’  entropy 
presented  in  Equations  (9)  and  (10)  is  mV  K-1.  Since  Ssc  is  depen¬ 
dent  on  the  electric  potential  (as  well  as  solid  phase  current  density 


By  setting  0s>p(L,t)  =  0.5exp(-0.07t)  V  as  a  boundary  condition  to 
the  porous  electrode  model,  the  simulated  short-circuit  current 
density  of  the  battery  is  shown  in  Fig.  6.  A  sudden  drop  of  the 
current  density  in  a  few  seconds  after  short-circuit  takes  place  is 
observed.  When  the  time  is  approximately  100  s,  the  sudden  drop 
bends  into  a  more  gradual  linear  declination,  then  continuing  in 
this  way  until  zero  is  reached  or  calculation  is  terminated. 

Fig.  7  shows  a  comparison  of  the  simulated  short-circuit  heat 
generation  rate  (Ssc)  when  0s>p(L,t)  =  0.5  V  and  when 
0s, p (L,t)  =  0.5exp(-0.07t)  V.  Although  the  overall  profiles  of  the  two 
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4.31274309  exp(0.571536523y)  +  1.281681122  sin(-4.9916739y)  -  0.090453431  sin(-20.9669665y  +  12.5788250) 
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4.14532933  +  8.1471 13434y  -  26.064581y2  +  12.7660158y3 


(10) 


and  SOC)  of  the  electrodes,  two  distinct  simulated  short-circuit 
heat  generation  rates  obtained  by  setting  two  different  types  of 
0s, p (L,t)  are  expected.  A  comparison  of  the  difference  is  shown  in 
Section  3.3. 

3.2.  Electrochemical  parameters  validation 

Various  electrochemical  parameters  need  to  be  determined 
in  order  to  apply  the  porous  electrode  theory  for  modeling 
battery  short-circuit.  Although  some  of  the  parameters  of  LMO 
cells  are  provided  by  Ref.  [22],  there  are  still  numerous  elec¬ 
trochemical  quantities  that  are  specific  to  each  cell,  such  as  the 
SEI  film  resistance  and  the  electrodes’  volume  fractions.  These 
values  could  be  obtained  by  comparing  the  simulated  galva- 
nostatic  discharge  curves  with  those  obtained  from  our 
experiment. 

In  this  experiment,  the  cell  was  discharged  inside  an  oven  of 
25.7  °C  ambient  temperature  by  connecting  the  battery  to  a  charge/ 
discharge  station  through  electric  wires  that  pass  through  the 
oven’s  door.  The  battery’s  output  voltage  and  current  during  the 
test  were  carefully  monitored  and  controlled  by  the  equipment, 
respectively,  to  obtain  reliable  galvanostatic  discharge  curves  of  the 
cell. 

On  the  other  hand,  the  simulated  discharge  curves,  recording 
the  evolution  of  0sp(x,  t)\.x=L  against  the  cell  capacity,  can  be 
plotted  by  solving  the  14  time-  and  space-dependent  variables  (i.e., 


curves  are  similar,  0s>p(L,t)  =  0.5  V  leads  to  an  underestimation  of 
the  maximum  value  by  19.6%  compared  with  that  obtained  by 
0S,P (L,t)  =  0.5exp(-0.07t)  V.  Therefore,  in  order  to  carry  out  a  more 
credible  short-circuit  simulation  model,  we  believe  that  a  precise 


Fig.  4.  The  simulated  (solid  lines)  and  experimental  (markers)  discharge  curves  of  the 
ITRI  5.25  Ah  soft  package  prismatic  cell.  The  effective  interfacial  area  of  cell,  Acen,  is 
taken  as  0.537  m2. 
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Table  1 

Electrochemical  parameters  for  simulating  the  ITRI  5.25  Ah  LMO  cell. 


Parameter 

Negative 

Electrolyte 

Positive 

electrode 

electrode 

Ds,i  (m2  s— 1 ) 

2  x  10  13a 

— 

1  x  10~13a 

Rf,,-  (Q  m2) 

0.06e 

— 

0 

vs,i  (S  m-1) 

100a 

— 

3.8a 

Rs,i  (fim) 

lle 

— 

8.5a 

,Z 

0.41 0e 

— 

0.297e 

£\,i 

0.357a 

1 

0.444a 

cM  (mol  m  3) 

17644m 

— 

11636m 

c° j-  (mol  ur3) 

13223m 

— 

228.1 6m 

kSj  (m2-5  mol-0-5  s"1) 

2  x  10  lle 

— 

2  x  10~lle 

U  (pm) 

74m 

20m 

107m 

Dj  (m2  s"1) 

3.227  x  10~10e 

n\  (S  ur1) 

0.914a 

t+ 

0.363 

Tfc  (°C) 

25 

eEstimated. 
mMeasured. 
a  Ref  [22]. 


Fig.  6.  Assumed  positive  electrode’s  electric  potential  0SjP(L,t)  =  0.5exp(-0.07t)  V 
(black  line)  and  the  corresponding  current  density  (red  line)  of  the  short-circuit  cell 
when  the  temperature  is  kept  constant  (25.7  °C)  in  the  simulation.  The  inset  is  an 
enlarged  plot  of  the  two  curves  around  the  initial  of  the  short  circuit.  (For  interpre¬ 
tation  of  the  references  to  color  in  this  figure  legend,  the  reader  is  referred  to  the  web 
version  of  this  article.) 


Fig.  7.  Comparison  of  the  simulated  short-circuit  heat  generation  rates  (Ssc)  by 
applying  the  two  different  positive  electrode’s  electric  potential  (0S>P)  at  25.7  °C.  The 
significant  difference  of  the  maximum  values  indicates  that  a  precise  description  of  0s  p 
is  critical. 


illustration  of  short-circuit  electric  potential  boundary  condition  is 
essential. 

Other  than  the  electric  potential  profile,  the  influences  of  SOC  on 
battery  heating  can  also  be  revealed  by  simulating  short-circuit 
with  the  porous  electrode  theory.  As  demonstrated  in  Fig.  8,  the 
maximum  value  of  Ssc  increases  while  SOC  decreases,  indicating 
that  a  more  severe  exothermic  reaction  occurs  with  smaller  SOCs. 
However,  the  reaction  decays  to  zero  more  rapidly  in  comparison  to 
a  fully  charged  cell  (SOC  =  1),  suggesting  that  the  time  interval  of 
exothermic  reaction  reduces  with  SOC. 


4.  Modeling  battery  thermal  runaway 

The  oven  test  is  frequently  conducted  to  observe  the  safety  of 
the  cell  under  high  temperature  exposure.  It  reveals  the  overall 
thermal  behavior  of  LIBs  under  thermal  abuse  condition.  To  give  a 
complete  illustration  of  the  battery’s  thermal  runaway  behavior 
when  it  is  exposed  to  elevated  temperature,  we  adopted  the  ther¬ 
mal  abuse  equations  17,18]  to  simulate  the  thermal  response  of  nail 
penetrated  LIBs  after  thermal  decomposing  reactions  in  the  cell  is 
triggered. 


Fig.  8.  Short-circuit  heat  generation  rates  (Ssc)  of  cells  short-circuited  at  different  SOCs 
when  T  =  25.7  °C.  Cells  with  smaller  SOCs  tend  to  have  greater  Ssc  peak  values,  but  they 
also  decay  more  rapidly  when  compared  with  cells  with  greater  SOCs.  The  inset  is  an 
enlargement  around  the  peaks  of  the  curves. 
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4.2.  Model  setup 


4.2.  Characteristics  of  heat  generation 


The  setup  of  the  thermal-runaway  model  considers  four  re¬ 
actions  of  the  jelly-roll  materials  observed  in  previous  studies:  the 
SEI  decomposition,  the  anode-electrolyte  reaction,  the  cathode- 
electrolyte  reaction,  and  the  electrolyte  decomposition  reaction 
[17,18  .  The  characteristics  of  the  reactions  are  dependent  on  their 
distinct  material  properties.  Therefore,  they  are  individually 
described  by  the  four  Arrhenius-type  equation  sets  shown  as 
follows: 

I  SEI  decomposition 


•^sei  —  ^sei^c^sei 


'a, sei 


RTi 


ta 


* 

C  • 


(ii) 


dcsei 

dt 


=  -Asei  exp 


'a, sei 


RT\ 


ta 


'ser 


(12) 


IIAnode-electrolyte  reaction 


Sn  —  HnWcAn  exp 


exp 


■£a,n 


RTi 


ta 


r 


(13) 


dc 


dt 


r  =  —An  exp 


-z 


exp 


— £a,n 


RT 


ta 


* 

r 


(14) 


dz 

dt=A"  exp 


-z 


exp 


— £a,n 


RT 


ta 


(15) 


III  Cathode-electrolyte  reaction 


Sp  —  HpWpApa(l  —  a:)  exp 


'3,P 


RT 


ta 


(16) 


da: 

dt 


Rp  =  Apa(l  -  a)exp 


'a,p 


RT 


ta 


(17) 


IV  Electrolyte  decomposition 


With  the  help  of  the  parameters  and  initial  values  shown  in 
Tables  2  and  3,  respectively,  the  four  heat  generation  rate  profiles  of 
the  thermal  abuse  equation  set,  i.e.  Ssei,  Sn,  Sp,  and  Se,  can  be  ob¬ 
tained.  In  order  to  clarify  the  different  temperatures  that  initiate 
the  exothermic  reactions,  we  give  a  temperature  sweep  of  1  °C  s_1 
to  observe  the  onset  of  the  action.  Simulation  results  of  the  ITRI 
5.25  Ah  LMO  cell  are  shown  in  Fig.  9,  and  the  heat  generation 
characteristics  of  the  four  reactions,  including  triggering  tempera¬ 
ture  and  the  whole  profile,  can  be  clearly  demonstrated.  In  addi¬ 
tion,  since  the  positive  electrode’s  material  is  LMO,  a  greater  peak  of 
Sp,  compared  to  [17,18],  is  obtained  based  on  the  material  param¬ 
eters  provided  by  Ref.  43  . 

From  the  result  of  assigning  temperature  sweep,  we  select  a 
specific  heat  generation  rate  of  105  W  m  3,  above  which  the  indi¬ 
vidual  reactions  are  triggered.  Due  to  its  unstable  nature,  the  SEI 
decomposition  occurs  at  the  lowest  temperature  value  of  the  four 
reactions  at  128  °C,  as  shown  in  Fig.  9.  Closely  following  behind  is 
the  anode-electrolyte  reaction  with  the  triggering  temperature  of 
146  °C.  The  cathode-electrolyte  reaction  takes  place  at  around 
182  °C,  showing  the  most  severe  exothermic  reaction,  and  the 
electrolyte  decomposition  occurs  last  at  approximately  229  °C.  It 
must  be  noted  that  different  LIBs  ought  to  have  different  heat 
generation  characteristics.  Batteries  with  various  shapes  and  sizes, 
or  even  with  novel  safety  designs,  are  also  critical  to  the  determi¬ 
nation  of  battery  thermal  runaway.  By  altering  the  frequency  factor 
Eaj  or  the  activation  energy  A/,  one  can  modulate  the  triggering 
temperature  and  the  overall  exothermic  reaction  profile  of  indi¬ 
vidual  batteries  to  fit  the  results  obtained  from  Differential  Scan¬ 
ning  Calorimetry  (DSC)  or  other  tests. 

5.  Modeling  nail  penetration 

Generally,  the  nail  penetration  process  of  LIBs,  which  involves  a 
dramatic  change  in  temperature,  is  the  coalescence  of  two  major 
mechanisms:  short-circuit  of  LIB,  where  it  happens  near  the  nail- 
battery  interface,  and  the  decomposition  of  jelly-roll  materials, 
occurring  within  the  whole  cell.  They  were  respectively  described 
in  Sections  3  and  4.  However,  to  provide  a  more  thorough  modeling 
framework  of  the  nail  penetration  test,  we  construct  a  compact 
model  that  simulates  the  nail  penetration  test  of  LIBs  by  coupling 
the  two  major  mechanisms. 


Se 


HeWeAe  exp 


— £a,e 


(18) 


dc* 

dt 


-Ae  exp 


— £a,e 

~RT^ 


* 


(19) 


The  readers  are  referred  to  [17,18]  for  a  detailed  description  of 
the  physical  interpretations  of  the  four  sets  of  reaction  kinetics. 


5.2.  Model  setup 

In  order  to  simulate  the  battery  temperature  distribution,  a 
three-dimensional  geometry  of  the  ITRI  5.25  Ah  soft  package 
prismatic  cell  nail  penetration  test  is  constructed,  and  the  sche¬ 
matic  diagram  is  shown  in  Fig.  10.  Several  simplifications  and  as¬ 
sumptions  of  the  original  complex  geometry  are  made  as  they  are 
peripheral  factors  or  do  not  have  significant  effect  to  the  modeling, 
as  stated  in  the  following.  The  current  collector  tabs  and  the  spirally 


Table  2 

Parameters  for  the  thermal  abuse  model  equations  obtained  from  Refs.  [17,18]. 


Parameter 

SEI 

decomposition 

Anode- 

electrolyte 

reaction 

Cathode- 

electrolyte 

reaction 

Electrolyte 

decomposition 

H.Ug-1) 

257 

1714 

400a 

155 

W,  (kg  m"3) 

6.104  x  102 

6.104  x  102 

1.438  x  103e 

4.069  x  102 

1.667  x  1015 

2.5  x  1013 

6.667  x  1013 

5.14  x  1025 

Ei  (J  mol'1) 

1.3508  x  105 

1.3508  x  105 

1.396  x  105 

2.74  x  105 

Estimated. 
a  Ref  [43]. 


Table  3 

Initial  values  for  thermal  abuse  model  variables 

[17,18]. 


Variable 

Initial  value 

C  sei 

0.15 

C  n 

0.75 

C*e 

1 

Z 

0.033 

a 

0.04 

Tea  (°C) 

25.7 
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Fig.  9.  The  heat  generation  rate  profiles  of  the  thermal  abuse  model,  i.e.  Ssei,  Sn,  5P,  and 
Se,  by  applying  a  temperature  sweep  of  1  °C  s_1.  The  cathode-electrolyte  reaction  (Sp) 
manifests  the  most  severe  exothermic  reaction  of  the  four.  The  inset  is  an  enlargement 
around  the  peaks  of  the  four  rates. 


Battery 


Fig.  10.  Schematic  diagram  of  the  three  dimensional  nail  penetration  model  setup. 


wound  geometry  of  the  jelly-roll  inside  the  cell  are  neglected.  The 
prismatic  cell  is  treated  as  a  bulk  and  the  thermal  properties  are 
volume-averaged  values.  Since  the  packaging  material  is  thin,  the 
aluminum  foil  wrapping  is  disregarded.  No  dilation  of  the  prismatic 
cell  is  assumed  due  to  the  two  fixed  Bakelite  in  the  actual  test  and 
the  difficulty  of  constructing  a  deformable  geometry.  The  influence 
of  the  steel  bolts  for  fastening  the  two  Bakelite  is  not  considered.  All 
the  materials  used  in  the  three-dimensional  geometry  are  regarded 
as  homogeneous  and  isotropic,  with  their  thermal  properties 
shown  in  Table  4. 

The  nail  diameter  is  a  geometry  feature  that  cannot  be  neglec¬ 
ted,  since  it  determines  the  total  area  of  short-circuit  in  the  nail 
penetration  test.  As  can  be  expected,  a  larger  short-circuit  area 
causes  a  greater  risk  of  inflammation  or  explosion.  Moreover, 
studies  have  pointed  out  that  the  nail  can  also  dissipate  consider¬ 
able  amount  of  heat,  which  may  retard  the  temperature  rising  rate 
of  the  cell  and  reduce  the  risk  of  thermal  runaway  [39,40].  When  a 
cylindrical  nail  is  penetrated  into  the  LIB,  short-circuit  is  assumed 
to  occur  only  at  the  interface  between  the  penetrated  nail  and  the 
battery.  Hence,  it  is  plausible  to  propose  the  short-circuit  heat 
generation  rate  per  unit  area,  S^frf,  as 

Sscrf  =  ^nailSsc,  (20) 


Table  4 

Thermal  properties  of  objects  used  in  the  nail  penetration  test. 


Parameter 

Nail 

Bakelite 

Jelly-roll 

kf  (W  m  1  K"1) 

44.5 

0.385 

0.034a 

Pi  (kg  m  3) 

7850 

1318 

1700a 

Cp.iOkg-1  K-') 

475 

14.4 

830a 

a  Ref  [17]. 


which  transforms  the  heat  generation  per  unit  volume  to  that  per 
unit  area  of  the  interface.  This  S^rf  will  serve  as  the  boundary  heat 
source  at  the  nail  and  battery  junction  in  our  modeling. 

With  the  heat  generation  rate  from  battery  short-circuit  and 
thermal  abuse  equations,  the  temperature  distribution  of  the  con¬ 
structed  three-dimensional  geometry  is  calculated  by  thermal  en¬ 
ergy  balance  between  the  generated  and  dissipated  heat  in  the  i-th 
mesh  element  [41,42] : 


ftCp,,g2  =  -v.(k,vr,)  +  y>,  (2i) 

j 

where  j  =  sei,  n,  p,  and  e,  and  V  is  the  gradient  operator  for  the 
three-dimensional  geometry’s  Cartesian  coordinate.  Furthermore, 
convective  cooling  is  applied  to  all  the  geometry’s  surfaces  that  are 
exposed  to  the  ambient: 


qconv  =  h  (lsurf  _  Tamb) ,  (22) 

where  h  is  the  heat  transfer  coefficient  and  the  ambient  tempera¬ 
ture  7*mb  is  assumed  to  remain  constant  at  the  initial  value  (25.7  °C) 
throughout  the  test. 

Through  Equation  (21),  the  equations  for  simulating  battery 
short-circuit  and  thermal  runaway  could  therefore  be  coupled 
together.  Their  connection  is  shown  in  Fig.  11,  where  short-circuit 
and  thermal  abuse  equations  calculate  their  corresponding  heat 
generation  rates,  and  the  time-space  variation  of  cell  temperature  is 
evaluated  via  thermal  energy  balance.  The  calculated  temperature 
is  treated  as  a  variable  that  is  substituted  back  into  the  corre¬ 
sponding  equations,  where  they  calculate  and  return  the  updated 
heat  generation  rates  to  Equations  (21)  and  (22),  allowing  the  cell 
temperature  distribution  of  the  next  time  step  to  be  evaluated. 

Although  the  simulation  of  battery  short-circuit  and  thermal 
abuse  require  the  input  of  temperature,  it  should  be  noted  that  the 
heat  generated  by  battery  short-circuit  is  applied  only  at  the  nail- 
battery  interface,  whereas  thermal  abuse  is  in  accordance  with 
the  battery’s  bulk  geometry.  Accordingly,  the  temperature  values 
for  the  two  equation  sets  should  have  distinct  locations  as  their 
source.  The  temperature  for  battery  short-circuit,  Tso  is  taken  as  the 
domain-averaged  value  at  the  nail-battery  interface.  As  will  be 
shown  later,  rsc  has  a  maximum  value  of  over  400  °C  during  the  nail 
penetration  process.  On  the  other  hand,  it  has  been  shown  that 
once  the  temperature  at  arbitrary  location  in  the  jelly-roll  exceeds 
the  material  decomposition  triggering  value,  a  whole  cell  thermal 
runaway  is  likely  to  occur  9].  As  a  result,  the  temperature  for 
thermal  abuse  model,  Tta,  is  taken  as  the  maximum  value  of  the 
whole  cell  to  probe  for  any  temperature  that  can  activate  the 
thermal  abuse  reactions. 


Fig.  11.  Coupling  between  the  short-circuit  and  the  thermal  abuse  modelings  via 
thermal  energy  balance  equation  in  the  nail  penetration  simulation. 
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5.2.  Simulation  result 

By  coupling  battery  short-circuit  and  thermal  abuse  through  the 
thermal  energy  balance  equation,  the  simulation  result  of  the  ITRI 
5.25  Ah  soft  package  prismatic  cell  nail  penetration  test  is  pre¬ 
sented  in  Fig.  12.  As  shown  in  the  diagram,  a  large  current  density 
from  the  battery  is  observed  due  to  short-circuit  occurrence  in  the 
test.  However,  the  current  drops  to  almost  20%  of  its  initial  value 
after  400  s,  because  of  the  depleting  active  material  initially 
intercalated  in  the  anode.  Coexisting  with  the  large  output  current 
is  the  severe  short-circuit  exothermic  reaction  at  the  nail-battery 
interface,  which  includes  the  Joule  heating  and  the  entropic  heat¬ 
ing  of  lithium  redox  reaction  as  shown  in  Equation  (8).  The  short- 
circuit  exothermic  heat  generation  rate  Ssc  reaches  its  maximum 
value  at  around  125  s  and  decays  gradually  with  time,  estimated  to 
also  result  from  depletion  of  the  active  material.  An  abrupt  drop  of 
current  density  and  jump  of  heat  generation  rate  at  50  s  are  also 
observed  in  the  figure.  The  cause  of  this  unexpected  variation  is 
explained  later  in  this  section  with  further  observation  of  other 
variables  in  the  model. 

During  the  initial  stage  with  only  the  occurrence  of  nail  pene¬ 
tration  and  no  thermal  runaway,  short-circuit  heat  generation  rate 
provides  the  initial  heating  of  the  cell  at  the  nail-battery  interface. 
Accordingly,  the  two  distinct  temperatures  for  the  two  models,  rsc 
and  Tta,  rise,  as  presented  in  Fig.  13.  The  two  temperatures  have  an 
overall  similar  behavior  within  the  interval  of  calculation.  However, 
Tta  instantly  jumps  to  approximately  450  °C  at  50  s  while  Tsc  con¬ 
tinues  in  a  more  gradual  path  to  its  maximum  value  of  just  above 
400  °C.  The  departure  of  the  two  temperatures  ends  at  about  75  s, 
and  continues  in  a  nearly  identical  path  until  the  end  of  calculation. 
The  instant  jump  of  the  temperature  at  50  s  could  be  explained  by 
observing  Fig.  14,  which  illustrates  the  thermal  runaway  heat 
generation  rate  profiles  with  specifying  Tta  as  its  temperature.  By 
comparing  this  result  with  the  one  that  was  not  coupled  (Fig.  9),  the 
thermal  runaway  heat  generation  rate  profiles  based  on  the  Tta 
value  (Fig.  14)  are  more  concentrated  than  those  of  temperature 
sweep  in  the  oven  test,  indicating  that  main  part  of  the  exothermic 
reaction  happens  within  a  short  time  interval,  specifically  at  around 
50  s.  This  characteristic  reveals  the  reason  for  the  unexpected 
variations  at  50  s  of  short-circuit  current  density  and  heat  gener¬ 
ation  rate  in  Fig.  12.  It  can  also  explain  the  departure  of  Tsc  and  Tta  at 
the  same  second,  since  it  is  the  moment  when  whole  cell  thermal 
runaway  occurs. 

Simulated  spatio-temporal  variation  of  cell  temperature  during 
the  nail  penetration  process  is  demonstrated  in  Fig.  15.  At  45  s,  the 
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Fig.  12.  Simulated  short-circuit  heat  generation  rate  (Ssc)  and  current  density  of  the 
nail  penetration  test  by  coupling  the  battery  short-circuit  and  thermal  abuse  equations. 
The  nail  penetration  occurs  at  0  s. 
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Fig.  13.  Comparison  of  Tsc  and  Tta.  The  two  curves  present  a  similar  profile  until  about 
50  s,  indicating  the  occurrence  of  an  additional  exothermic  reaction  other  than  short 
circuit. 


heating  of  the  cell  is  contributed  by  short-circuit  heat  generation 
rate,  forming  a  local  hot  spot  at  the  point  of  penetration.  Thermal 
runaway  occurs  at  50  s,  causing  the  temperature  of  the  whole  cell  to 
rise  abruptly.  Due  to  the  heat  conduction  lag,  maximum  surface 
temperature  of  the  cell  is  detected  at  60  s.  After  60  s,  convection 
mechanism  determines  the  temperature  variation  of  the  cell  since 
the  thermal  abuse  heat  generation  rate  dies  out  quickly  after  its 
short  occurrence. 

Referring  to  the  experimental  results,  we  capture  the  simulated 
temperature  time  variation  data  at  spots  corresponding  to  the 
experiment  thermal  couple  locations.  The  compared  simulation 
and  experiment  results  are  shown  in  Fig.  16.  Since  thermal  couples 
A  and  B  (as  well  as  C  and  D)  are  located  at  symmetric  positions,  we 
only  present  the  data  obtained  from  one  of  the  spots  to  represent 
our  simulation  result.  To  verify  the  obtained  temperature  profiles 
from  simulation,  there  are  three  indexes  that  can  be  compared 
with  the  experimental  data.  First,  the  slopes  of  the  temperatures  at 
thermal  runaway  obtained  from  the  experiments  are  nearly  ver¬ 
tical.  The  simulation  traces  out  this  characteristic  successfully 
since  most  of  the  heat  generation  reaction  concentrates  in  a  nar¬ 
row  time  interval  (Fig.  14).  Second,  temperature  from  thermal 
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Fig.  14.  Simulated  heat  generation  rate  profiles  of  the  thermal  abuse  model  by  inte¬ 
grating  short-circuit  electrochemical  model  and  the  thermal  abuse  model.  Almost  all 
the  exothermic  reactions  are  concentrated  in  a  short  time  interval,  pointing  out  the 
occurrence  of  a  severe  whole  cell  thermal  runaway  at  the  instant.  The  inset  is  the 
enlarged  plot  around  the  moment  of  thermal  runaway. 
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Fig.  15.  Simulated  spatio-temporal  temperature  in  the  nail  penetration  test. 


couple  E,  which  is  nearest  to  the  penetration  point,  peaks  at 
310.6  °C.  The  simulation  reveals  a  maximum  temperature  of 
327.3  °C,  a  5.38%  error,  from  the  temperature  profile  corresponding 
to  thermal  couple  E.  Third,  the  profiles  of  temperature  drop  after 
the  apex  are  primarily  determined  by  the  ambient  conditions. 
Since  ventilation  of  the  chamber  is  turned  on  after  the  penetration, 
a  convection  cooling  will  occur.  A  cooling  boundary  condition, 
which  is  then  set  with  the  heat  transfer  coefficient  at  a  higher 
value  of  h  =  71.7  W  rrT2  K-1  (compared  with  h  ~  7.6  W  m~2  K_1 
for  natural  convection),  produces  a  temperature  profile  that  fits 
well  with  the  experimental  results. 

6.  Conclusion 

We  have  proposed  a  novel  method  by  using  the  porous  elec¬ 
trode  theory  to  model  the  electrochemical  processes  of  short- 
circuited  LIBs.  By  doing  so,  we  have  calculated  the  output  current 
of  short-circuit  batteries  and  predicted  the  heat  generation  rates  of 
cells  with  different  electrochemical  parameters,  such  as  SOC.  In 
addition,  with  the  battery  thermal  abuse  equations  from  the  oven 
test,  we  have  also  simulated  the  thermal  runaway  behavior  of  LIBs. 


Fig.  16.  Comparison  of  the  (a)  simulated  and  (b)  experimental  results.  The  results 
show  good  agreement  on  describing  the  temperature  rise,  the  maximum  temperature, 
and  the  temperature  declination  due  to  ventilation. 


By  combining  the  two  modelings,  our  simulation  result  of  cell 
temperature  variation  shows  great  similarity  with  the  nail  pene¬ 
tration  test.  Therefore,  we  believe  this  method  is  capable  to  accu¬ 
rately  model  LIB  nail  penetration  tests.  Furthermore,  our  model  can 
be  extended  to  predict  the  effectiveness  of  LIBs  with  certain  safety 
features,  such  as  cells  with  additives  that  can  shut  down  the  mass 
or  charge  transfer  process  when  the  battery  is  short-circuited. 
Finally,  regarding  cases  where  a  nail  is  penetrated  into  cells  with 
new  electrodes/electrolytes  or  cells  that  are  cycled,  preliminary 
experiments  are  required  to  obtain  the  parameters  of  cells  with 
new  materials,  and  additional  models  need  to  be  constructed  to 
describe  jelly-roll  deterioration  after  cycling  tests.  Future  work  will 
explore  these  issues  in  more  detail. 
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k:  heat  transfer  coefficient,  W  m-1  K-1 

ks:  reaction  rate  constant,  m25  mol-0-5  s-1 

L:  thickness,  m 

rnau:  nail  radius,  m 

R:  gas  constant,  8.314  J  K-1  mol-1 

Rf.  SEI  film  resistance,  Q  m2 

Rs:  radius  of  electrode  particle,  pm 

S:  heat  generation  rate,  W  m-3 

T:  temperature,  K 

t:  time,  s 

t+:  transference  number 

Uc:  open-circuit  potential  at  25  °C,  V 

W:  active  content  density,  kg  m-3 

x:  state  of  charge  of  negative  electrode 

y:  state  of  charge  of  positive  electrode 

z;  dimensionless  SEI  thickness 

a:  degree  of  conversion 

ef  Filler  volume  fraction 

ep  electrolyte  volume  fraction 

es:  porous  electrode  volume  fraction 

<fp  electric  potential  of  electrolyte  phase,  V 

(j)s:  electric  potential  of  solid  phase,  V 

y;  Bruggeman  exponential  coefficient,  1.5 

r\s:  surface  over-potential,  V 

up  electrolyte  phase  ionic  conductivity,  S  m 

p:  density,  kg  m-3 

as:  solid  phase  electric  conductivity,  S  m-1 


Nomenclature 


Subscript 


a:  specific  interfacial  area,  m2  m-3 

A:  frequency  factor,  s_1 

Aceii:  total  effective  interfacial  area  of  cell,  m2 

c:  concentration  of  salt  in  electrolyte,  mol  m-3 

cs:  concentration  of  lithium  in  electrode,  mol  m-3 

ce:  dimensionless  concentration  of  electrolyte 

c*:  dimensionless  amount  of  Li  intercalated  within  the  carbon 

cse[:  dimensionless  amount  of  SEI  Li-containing  metastable  species 

Cp:  heat  capacity,  J  kg-1  K_1 

Dp  electrolyte  phase  lithium-ion  diffusivity,  m2  s-1 

Ds:  solid  phase  lithium  diffusivity,  m2  s-1 

Ea:  activation  energy,  J  mol-1 

F:  Faraday  constant,  96485  C  mol-1 

H:  specific  energy,  J  g-1 

h:  heat  transfer  coefficient,  W  m-2  K-1 

Ip  electrolyte  phase  current  density,  A  m-2 

Is:  solid  phase  current  density,  A  m-2 

jiP  pore  wall  flux  of  lithium-ion,  mol  m-2  s-1 


e:  separator  region 

i:  for  component,  species,  or  reaction  i 

n:  negative  electrode 

p;  positive  electrode 

sc:  short-circuit 

ta:  thermal  abuse 

cell:  battery 

sei:  solid-electrolyte  interface 
Superscript 


0:  initial  state 
M:  maximum 
eff:  effective 
amb:  ambient 
surf:  surface 


